home *** CD-ROM | disk | FTP | other *** search
/ The Datafile PD-CD 1 Issue 2 / PDCD-1 - Issue 02.iso / _utilities / utilities / 001 / meschach / !Meschach / h / iter < prev    next >
Text File  |  1994-03-08  |  7KB  |  249 lines
  1.  
  2. /**************************************************************************
  3. **
  4. ** Copyright (C) 1993 David E. Steward & Zbigniew Leyk, all rights reserved.
  5. **
  6. **                 Meschach Library
  7. ** 
  8. ** This Meschach Library is provided "as is" without any express 
  9. ** or implied warranty of any kind with respect to this software. 
  10. ** In particular the authors shall not be liable for any direct, 
  11. ** indirect, special, incidental or consequential damages arising 
  12. ** in any way from use of the software.
  13. ** 
  14. ** Everyone is granted permission to copy, modify and redistribute this
  15. ** Meschach Library, provided:
  16. **  1.  All copies contain this copyright notice.
  17. **  2.  All modified copies shall carry a notice stating who
  18. **      made the last modification and the date of such modification.
  19. **  3.  No charge is made for this software or works derived from it.  
  20. **      This clause shall not be construed as constraining other software
  21. **      distributed on the same medium as this software, nor is a
  22. **      distribution fee considered a charge.
  23. **
  24. ***************************************************************************/
  25.  
  26.  
  27. /* iter.h  14/09/93 */
  28.  
  29. /* 
  30.  
  31.   Structures for iterative methods
  32.  
  33. */
  34.  
  35. #ifndef ITERHH
  36.  
  37. #define ITERHH
  38.  
  39. /* RCS id: $Id: iter.h,v 1.2 1994/03/08 05:48:27 des Exp $  */
  40.  
  41.  
  42. #include    "sparse.h"
  43.  
  44.  
  45. /* basic structure for iterative methods */
  46.  
  47. /* type Fun_Ax for functions to get y = A*x */
  48. #ifdef ANSI_C
  49. typedef VEC  *(*Fun_Ax)(void *,VEC *,VEC *);
  50. #else
  51. typedef VEC *(*Fun_Ax)();
  52. #endif
  53.  
  54.  
  55. /* type ITER */
  56. typedef struct Iter_data {
  57.    int shared_x;   /* if TRUE then x is shared and it will not be free'd */ 
  58.    int shared_b;   /* if TRUE then b is shared and it will not be free'd */
  59.    unsigned k;   /* no. of direction (search) vectors; =0 - none */
  60.    int limit;    /* upper bound on the no. of iter. steps */
  61.    int steps;    /* no. of iter. steps done */
  62.    Real eps;     /* accuracy required */
  63.    
  64.    VEC *x;       /* input: initial guess;
  65.             output: approximate solution */
  66.    VEC *b;       /* right hand side of the equation A*x = b */
  67.  
  68.    Fun_Ax   Ax;         /* function computing y = A*x */
  69.    void *A_par;         /* parameters for Ax */
  70.  
  71.    Fun_Ax  ATx;         /* function  computing y = A^T*x;
  72.                            T = transpose */
  73.    void *AT_par;         /* parameters for ATx */
  74.  
  75.    Fun_Ax  Bx; /* function computing y = B*x; B - preconditioner */
  76.    void *B_par;         /* parameters for Bx */
  77.  
  78. #ifdef ANSI_C
  79.  
  80. #ifdef PROTOTYPES_IN_STRUCT
  81.    void (*info)(struct Iter_data *, double, VEC *,VEC *);
  82.             /* function giving some information for a user;
  83.            nres - a norm of a residual res */
  84.    
  85.    int (*stop_crit)(struct Iter_data *, double, VEC *,VEC *);
  86.            /* stopping criterion:
  87.           nres - a norm of res;
  88.           res - residual;
  89.         if returned value == TRUE then stop;
  90.         if returned value == FALSE then continue; */
  91. #else
  92.    void (*info)();
  93.    int  (*stop_crit)();
  94. #endif /* PROTOTYPES_IN_STRUCT */
  95.  
  96. #else
  97.  
  98.    void (*info)();
  99.             /* function giving some information for a user */
  100.    
  101.    int (*stop_crit)();
  102.            /* stopping criterion:
  103.         if returned value == TRUE then stop;
  104.         if returned value == FALSE then continue; */
  105.  
  106. #endif /* ANSI_C */
  107.  
  108.    Real init_res;   /* the norm of the initial residual */
  109.  
  110. }  ITER;
  111.  
  112.  
  113. #define INULL   (ITER *)NULL
  114.  
  115. /* type Fun_info */
  116. #ifdef ANSI_C
  117. typedef void (*Fun_info)(ITER *, double, VEC *,VEC *);
  118. #else
  119. typedef void (*Fun_info)();
  120. #endif
  121.  
  122. /* type Fun_stp_crt */
  123. #ifdef ANSI_C
  124. typedef int (*Fun_stp_crt)(ITER *, double, VEC *,VEC *);
  125. #else
  126. typedef int (*Fun_stp_crt)();
  127. #endif
  128.  
  129.  
  130.  
  131. /* macros */
  132. /* default values */
  133.  
  134. #define ITER_LIMIT_DEF  1000
  135. #define ITER_EPS_DEF    1e-6
  136.  
  137. /* other macros */
  138.  
  139. /* set ip->Ax=fun and ip->A_par=fun_par */
  140. #define iter_Ax(ip,fun,fun_par) \
  141.   (ip->Ax=(Fun_Ax)(fun),ip->A_par=(void *)(fun_par),0)
  142. #define iter_ATx(ip,fun,fun_par) \
  143.   (ip->ATx=(Fun_Ax)(fun),ip->AT_par=(void *)(fun_par),0)
  144. #define iter_Bx(ip,fun,fun_par) \
  145.   (ip->Bx=(Fun_Ax)(fun),ip->B_par=(void *)(fun_par),0)
  146.  
  147. /* save free macro */
  148. #define ITER_FREE(ip)  (iter_free(ip), (ip)=(ITER *)NULL)
  149.  
  150.  
  151. /* prototypes from iter0.c */
  152.  
  153. #ifdef ANSI_C
  154. /* standard information */
  155. void iter_std_info(ITER *ip,double nres,VEC *res,VEC *Bres);
  156. /* standard stopping criterion */
  157. int iter_std_stop_crit(ITER *ip, double nres, VEC *res,VEC *Bres);
  158.  
  159. /* get, resize and free ITER variable */
  160. ITER *iter_get(int lenb, int lenx);
  161. ITER *iter_resize(ITER *ip,int lenb,int lenx);
  162. int iter_free(ITER *ip);
  163.  
  164. void iter_dump(FILE *fp,ITER *ip);
  165.  
  166. /* copy ip1 to ip2 copying also elements of x and b */
  167. ITER *iter_copy(ITER *ip1, ITER *ip2);
  168. /* copy ip1 to ip2 without copying elements of x and b */
  169. ITER *iter_copy2(ITER *ip1,ITER *ip2);
  170.  
  171. /* functions for generating sparse matrices with random elements */
  172. SPMAT    *iter_gen_sym(int n, int nrow);
  173. SPMAT    *iter_gen_nonsym(int m,int n,int nrow,double diag);
  174. SPMAT    *iter_gen_nonsym_posdef(int n,int nrow);
  175.  
  176. #else
  177.  
  178. void iter_std_info();
  179. int iter_std_stop_crit();
  180. ITER *iter_get();
  181. int iter_free();
  182. ITER *iter_resize();
  183. void iter_dump();
  184. ITER *iter_copy();
  185. ITER *iter_copy2();
  186. SPMAT    *iter_gen_sym();
  187. SPMAT    *iter_gen_nonsym();
  188. SPMAT    *iter_gen_nonsym_posdef();
  189.  
  190. #endif
  191.  
  192. /* prototypes from iter.c */
  193.  
  194. /* different iterative procedures */
  195. #ifdef ANSI_C
  196. VEC  *iter_cg(ITER *ip);
  197. VEC  *iter_cg1(ITER *ip);
  198. VEC  *iter_spcg(SPMAT *A,SPMAT *LLT,VEC *b,double eps,VEC *x,int limit,
  199.         int *steps);
  200. VEC  *iter_cgs(ITER *ip,VEC *r0);
  201. VEC  *iter_spcgs(SPMAT *A,SPMAT *B,VEC *b,VEC *r0,double eps,VEC *x,
  202.          int limit, int *steps);
  203. VEC  *iter_lsqr(ITER *ip);
  204. VEC  *iter_splsqr(SPMAT *A,VEC *b,double tol,VEC *x,
  205.           int limit,int *steps);
  206. VEC  *iter_gmres(ITER *ip);
  207. VEC  *iter_spgmres(SPMAT *A,SPMAT *B,VEC *b,double tol,VEC *x,int k,
  208.            int limit, int *steps);
  209. MAT  *iter_arnoldi_iref(ITER *ip,Real *h,MAT *Q,MAT *H);
  210. MAT  *iter_arnoldi(ITER *ip,Real *h,MAT *Q,MAT *H);
  211. MAT  *iter_sparnoldi(SPMAT *A,VEC *x0,int k,Real *h,MAT *Q,MAT *H);
  212. VEC  *iter_mgcr(ITER *ip);
  213. VEC  *iter_spmgcr(SPMAT *A,SPMAT *B,VEC *b,double tol,VEC *x,int k,
  214.           int limit, int *steps);
  215. void    iter_lanczos(ITER *ip,VEC *a,VEC *b,Real *beta2,MAT *Q);
  216. void    iter_splanczos(SPMAT *A,int m,VEC *x0,VEC *a,VEC *b,Real *beta2,
  217.                MAT *Q);
  218. VEC  *iter_lanczos2(ITER *ip,VEC *evals,VEC *err_est);
  219. VEC  *iter_splanczos2(SPMAT *A,int m,VEC *x0,VEC *evals,VEC *err_est);
  220. VEC  *iter_cgne(ITER *ip);
  221. VEC  *iter_spcgne(SPMAT *A,SPMAT *B,VEC *b,double eps,VEC *x,
  222.           int limit,int *steps);
  223. #else
  224. VEC  *iter_cg();
  225. VEC  *iter_cg1();
  226. VEC  *iter_spcg();
  227. VEC  *iter_cgs();
  228. VEC  *iter_spcgs();
  229. VEC  *iter_lsqr();
  230. VEC  *iter_splsqr();
  231. VEC  *iter_gmres();
  232. VEC  *iter_spgmres();
  233. MAT  *iter_arnoldi_iref();
  234. MAT  *iter_arnoldi();
  235. MAT  *iter_sparnoldi();
  236. VEC  *iter_mgcr();
  237. VEC  *iter_spmgcr();
  238. void  iter_lanczos();
  239. void  iter_splanczos();
  240. VEC  *iter_lanczos2();
  241. VEC  *iter_splanczos2();
  242. VEC  *iter_cgne();
  243. VEC  *iter_spcgne();
  244.  
  245. #endif
  246.  
  247.  
  248. #endif  /* ITERHH */
  249.